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In order to support the anytime abort requirements of a manned mission to the Moon, 
the vehicle abort capabilities for the translunar and circumlunar phases of the mission 
must be studied. Depending on the location of the abort maneuver, the maximum return 
time to Earth and the available propellant, two different kinds of return trajectories can 
be calculated: direct and fly-by. This paper presents a new method to compute these 
return trajectories in a deterministic and fast way without using numerical optimizers. 
Since no simplifications of the gravity model are required, the resulting trajectories are 
very accurate and can be used for both mission design and operations. This technique has 
been extensively used to evaluate the abort capabilities of the Orion/Altair vehicles in the 
Constellation program for the translunar phase of the mission. 


I. Introduction 

A bort return trajectories play an important role In the performance of a vehicle designed for manned 
missions. For a given nominal mission, different abort scenarios can be presented. In this paper, we will 
cover those scenarios related to the on-orbit performance part and only for Earth to Moon transfers. The 
main goal will be to calculate abort return trajectories, given the vehicle state and some Entry Interface (El) 
constraints, in a fast and accurate way. In order to evaluate the performance of the vehicle, abort trajectories 
will have to be calculated along the nominal trajectory at different times. In Figure 1, we can see two different 
abort scenarios for a typical manned mission to the Moon: (a) the abort maneuver is implemented before 
the Lunar-Orbit-Insertion maneuver (LOI) and (b) the abort maneuver is implemented after a failed LOI. 
In both cases two strategies can be implemented: a direct return and a Moon flyby. Choosing a direct or 
flyby return will depend on the abort scenario, distance to the Moon and operational constraints. 

In this paper, a method to compute high-fidelity feasible solutions to the abort return problem is de- 
scribed. The first part of the paper summarizes the two techniques used to compute the return trajectories: 
Pseudostate theory developed by S. Wilson 5,7 and the optimal El targeter developed by this author. 4 The 
second part discusses how to integrate these two techniques to calculate abort return trajectories in the 
translunar and circumlunar phases for the direct and fly-by cases. A discussion about the optimality of the 
obtained return trajectories, when different entry interface constraints are present, is also included. Finally, 
some numerical examples and comparisons to the results obtained by numerical optimization are described. 

II. Pseudostate theory review 

In this section, a brief description of the main results and algorithms of Pseudostate Theory (PT) for 
Earth-Moon transfers is presented. For a complete explanation of PT and its applications see Refs. 5 ’'. In 
general, in PT a three-body trajectory is approximated by two osculating conics at the beginning and at the 
end of the propagation (see Figure 2). For the abort problem we are trying to solve, we have an initial state 
around the Moon and a final state at the Earth El: 

initial state (Moon-centered inertial coordinate frame): f/,r/,v/ 
final state (Earth-centered inertial coordinate frame): tEi ,Re/ , V ei 
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Moon-centered inertial view 



(b) 


Figure 1. Abort trajectory scenario examples for cases: (a) pre-lunar-orbit-insertion abort and (b) post-lunar- 
orbit-insertion abort. Nominal Earth to Moon trajectory (red), abort return trajectories (black) and Earth 
trajectory (blue). 
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pseudostate sphere (rg) 



Figure 2. Pseudostate propagation for Earth-Moon transfers. Actual trajectory (red). Primary and secondary 
conics (solid black). 


The goal is to obtain an approximation of the El state resulting from the propagation of the initial state 
under the actual gravitational field (function Propagate): 

(R E i,V EI ) = Propagate (r /; v 7 , t 7 , t EI ) 

The two osculating conics at the above points are: 


secondary conic: (ts, rs, vg) = Conic2distance (r i, Vi, rs, Hi) (1) 

primary conic: (R E i, V ei) = Conic (R 7 *, VJ*, t E i - ti, he) (2) 

The secondary conic propagates the initial Moon-centered inertial state to a distance rs (function 
Conic2distance) using the gravitational constant hi- This distance rs is defined by the pseudostate sphere 
radius. For an Earth-Moon transfer rs is about 24 times the radius of the Earth. From the propagation of 
the secondary conic we obtain (fs, rs,vs), the time and state where the trajectory pierces the pseudostate 
sphere. The primary conic propagates the primary Earth-centered pseudostate (f 7 , R}*, VJ*) up to the 
EI point (function Conic) using the gravitational constant he- These two conics are related through two 
pseudostates: secondary (rj, vj) and primary (R|*, VJ*) by 


* 

v / 

= v s 

( 3 ) 

* 

r / 

= r s + (t 7 - t s ) v s 

( 4 ) 

R *j* 

= R-m + r 7 + <^ 7 

( 5 ) 

v** 

= V m + vj + 6j 

( 6 ) 


where: 


(Rm, V m ): 


Moon position and velocity at tj in Earth-centered inertial coordinates 
calibration vectors. To take into account high-order gravitational fields and 
third-body perturbations (e.g. from the SUN). 

The propagation process is as follows: 


1. Use Eq. 1 to obtain (ts,rs,vs). 
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2. Use Eqs. 3, 4, 5 and 6 to obtain the two pseudostates: (rj,vj) and (RJ*, V}*) 

3. Use Eq. 2 to obtain (R ei^ei) 

>)c * * 

Note that for propagation-only purposes, no calibration is performed and therefore: 6j = dj = 0. 

A. Three-body Lambert problem 

We can modify the propagation algorithm described in the previous section to define a three-body targeting 
algorithm. Given an initial position and time relative to the secondary body (tg,rg) and a final postion 
relative to the primary body ( t E , Rb), we need to calculate the initial and final velocities (vg, Vg). We also 
need to know the gravitational constants associated with the primary /tg; and secondary /xg bodies together 
with the size of the pseudostate sphere associated with the secondary body rg. In Algorithm 1, we can see 
a pseudocode version of the proposed procedure. Initially, the secondary pseudostate position and the two 
calibration vectors are set to zero (line 1 in the algorithm). After that, the pseudostate loop is executed (lines 
2 to 15). This loop can be divided into two parts: pseudostate propagation (lines 3 to 9) and calibration 
(lines 10 to 14). 

In line 3, the current estimate of the primary pseudostate position vector is computed from the Moon 
position at the initial time, the secondary pseudostate and the calibration position vector. In line 4, the 
primary conic is calculated by solving the Lambert problem from the primary pseudostate position vector 
R to the final Earth position vector Rg; for a transfer time tE — ti- In line 5, the velocity of the secondary 
pseudostate is calculated vj. From Eq. (3), this is also the velocity at the end of the secondary conic. In 
line 6, the secondary conic will be calculated. For that, the INRFV problem will be solved (see section B). 
Once the secondary conic is calculated, we can now get the position vector of the secondary pseudostate rj 
from Eq. (4). Finally, we will iterate until the value of rj has converged with an specified tolerance (line 9). 
Once the pseudostate propagation part has converged, the calibration part of the algorithm will be executed 
(lines 10 to 14). In line 10, a numerical propagation of the initial state: f/,17, vg until time tE is performed. 
This numerical propagation will incorporate high-order gravity models, third-body perturbations, etc. In 
general, it will be performed by one or several calls to ODE solvers. If after numerically propagating the 
solution, the final state I? 7 is close to the target Rg; within some tolerance, the solution has been found 
and the algorithm exits (line 11). Otherwise, the primary pseudostate vectors RJ*, VJ* will be obtained by 
propagating a conic backward in time from the numerically propagated final state H E ,V £ (line 12). Using 
Eqs. (5) and (6), we can now calculate the calibration vectors Sj and 5j (lines 13 and 14). 

B. The INRFV problem 

The problem we are trying to solve can be defined as (see Figure 3): given an initial position vector ri 
and final velocity V2 at some distance r2, calculate the initial velocity vi (for a complete description of the 
solution to this problem see Ref. 6 ). For the abort problem we are solving in this paper, we will assume that 
v 2 = rg (pseudostate sphere radius) and r 2 > rq. From Figure 3, we can define the angle between r 1 and v 2 
as: 

ip = 9 + /3 2 
cos ip = r J v 2 

sin ip = ± \Jl — cos 2 ip (7) 

This problem requires the solution of a quatic equation: 

A 2 x 4 + 2 ABx 3 + (B 2 + C- 2 A) x 2 - 2 Bx + 1 - C = 0 (8) 

where: 

A= r - 2 ^ 2 , B=( l-™£)smip, C = cos 2 ip, x = sinp 2 

M r i V m / 

We can solve this quartic equation analytically. For this work the routine PA05AD from HSL 2 has been 
used. The formulation of Eq. (8) introduces extraneous solutions, in order to eliminate them we have used 
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Algorithm 1 Pseudostate three-body Lambert problem 
Given: tj, t E , r j, R s , rs, He and /i/ 

Obtain: v/ and 

Initialization 

1: r j = 5j = ij = 0 


Pseudostate loot 
2: do 


3: RJ* = R m +r}+S*j 

4: (VJ*, Vg) = Lambert (RJ*, Rg, t E - G, He) 

5: v* = VJ* - V m - 6* 

6: (v /} r s , t s ) = INRFV (r,, vj, r s , Hi) 

7: rj' = rj 

8: rj = r s + (t/ - ts) vj 

9: if r J — r J > tolerance cycle 

Calibration 

10: (r £ , V E ) = Propagate (rj, vj, f/, tg) 

11: if Re - R B < tolerance exit 

12: (RJ*, V**) = Conic (R b , V E , t, - t E , He ) 

13: = RJ* - R,„ - rj 

14: <T, = VJ* - V m - vj 

15: end do 


the following criteria (see Ref. 5 6 for a complete explanation): If cosip > 0 r]"v 2 > 0, sin/3 2 is the smaller 

of the real positive solutions. If cos '! p < 0 rjv 2 < 0, sin/3 2 is the larger of the real positive solutions. 
Once sin f3 2 has been found we can compute vi and r 2 as follows: 


fi x v 2 

c = — T~ 

sin ip 

9 = atan2(smip,cosip) — sin/3 2 

r 2 = r 2 [cos 0r i + sin# (c x fi)] 


a -2^ 

+ 2 — 

[r 2 v 2 sm /? 2 ) 

P = 

r 2 

n 

H 

r 2 (1 

— cosd) 

r 2 r i sin# 


V 

to 

II 


( 9 ) 


r 2 - fr i 

Vi = 

9 

One final note about the problem. The sign of sin ip can be chosen arbitrarily (see Eq. (7)). Therefore 
there are two antiparallel choices for the direction of the angular momentum (see Eq. (9)). In one case r x 
will be located after periapse (see Figure 3 solid red line) and in the second case iq will be located before 
periapse (see Figure 3 dashed red line). This is important for the abort problem we are trying to solve 
because just by changing the direction of the angular momentum we can specify the abort solutions to do a 
Moon flyby or not (direct return). In Figure 1, we can see that in both scenarios: cases (a) and (b), we can 
always find a direct and a Moon flyby return. 
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Figure 3. INRFV problem definition. Depending on the sign of sin two different solutions can be obtained 
(red solid and red dashed curves). 

III. Three-body abort targeting problem 

In this section, the three-body problem formulated in Algorithm 1 will be modified so it can be used as a 
targeting algorithm for abort problems. The main difference between the proposed abort-targeting algorithm 
and Algorithm 1 is the use of a different method to calculate the primary conic. Instead of using a Lambert 
arc for this algorithm, an optimal two-body abort targeting algorithm developed by the author 4 will be used. 
In general, this algorithm will solve the following problem: Given an initial state (£, R, V), a transfer time 
At to Entry Interface (El) and a gravitational constant //,, calculate an optimal impulsive maneuver AV 
such that at El the following constraints are satisfied: 


equality constraints: tei, Jei 
inequality constraints: A ei, Lei , Azei 


El radius and flight-path-angle. 

El longitude, latitude and azimuth. 


As described in Ref. 4 , this targeting algorithm has some properties that makes it suitable for the problem 
we are trying to solve in this paper. In general, the algorithm is non-iterative (an iteration procedure is used 
when the El velocity with respect to the At is calculated but this is done only once) and it provides a global 
optimal solution and feasibility information about the problem. 

The three-body abort targeting algorithm (see Algorithm 2) is very similar to the one in Algorithm 1. 
Therefore, only the differences between them will be described. Whereas the Lambert routine only needs 
the position vectors at the beginning and at the end of the arc, the optimal abort algorithm described in 
this section will need a complete state at the initial time and some constraints at EI. In order to compute 
the full primary pseudostate (£/, RJ*, V))*) , the secondary pseudostate will have to be defined first. In line 
2 of Algorithm 2, the secondary pseudostate is initialized and in line 6, the primary pseudostate velocity is 
calculated. In lines 7 and 8, the primary conic is calculated such that it satisfies the EI constraits specified: 
tei, 'Jei at the least or tei, 7 ei and some inequality constraints on Xei, Lei, Azei- The computation of 
the primary conic and the calibration section of the algorithm (from line 9 to 19) are identical to the ones 
in Algorithm 1. One final difference between both targeting algorithms is that as Algorithm 2 converges 


0, whereas previously only 


to the solution AVj * 
initialization of the secondary pseudostate 


0 and 


Re/ — R 


EI 


0. Finally, the 


/! v / 


in lines 1 and 2 is not unique. Through numerical tests 
it was shown that rj = 0, vj = v/ performs well in the type of problems solved in this paper. It is also 
worth mentioning that as suggested in Ref. 5 , the use of a pre-periapse pseudostate instead of using always 
the post-periapse one (as in Algorithm 2) can reduce the number of iterations. In numerical results we have 
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Algorithm 2 Three-body abort targeting problem 
Given: tj, t E i , rj, v/, r EI , Jei, rs, He and hi 
Obtain: A vj and El conditions 

Initialization 

1: rj = 0 
2: vj = v/ 

3: g*=S*=0 


Pseudostate loot 
4: do 


5: 

6 : 

7: 

8 : 

9: 

10 : 

11 

12 : 

13 


R T = Rm +r*j+S*j 
VJ* = V m +vJ +5j 

(AV”, Rb/, Vbj) = Abort (RJ*, V}*, r E i, Jei, t E i — ti, he) 

y** ^ y** _|_ ^y** 

vj = VJ* - V m - 5* 

(v/, r s , ts) = INRFV (rj, vj, r s , p/) 

r *' _ r * 

T I ~ T I 

r} =r s + (tj - ts) vj 

if IrJ -rj'l > tolerance cycle 


Calibration 

14: (tl' ei , V EJ ^) = Propagate (r 7 , vj, t/, t EI ) 

15: if Rs/ — R S j < tolerance exit 

16: (RJ*, VJ*) = Conic ^R^;/ , V 1/ — t E i , 

17: 5J = RJ* - R,„ - rj 

18: <5/ = VJ* - V m - vj 


19: end do 

20: A v 7 = Vj — v 7 


seen that the pre-periapse pseudostate reduces the number of iterations in cases where the trajectory is just 
entering the pseudostate sphere. Once the trajectory is inside the pseudostate sphere both pseudostates (pre 
and post-periapse) perform similarly for the cases evaluated in this work. 

IV. Optimality of the solutions 

In this section, we will discuss the optimality of the solutions for the three-body and two-body targeting 
algorithms described above. The three-body version of the algorithm, which includes calls to the two-body 
algorithm, will be used when the pseudostate conditions apply. That is when the vehicle is inside the 
pseudostate sphere. Otherwise, the two-body version will be used. One difference in the performance of 
both algorithms is that the three-body version will generate high-fidelity solutions. Since the three-body 
version incorporates a calibration phase where actual trajectories are calculated, the final answer is always 
a return trajectory numerically integrated with high-order gravity fields, third-body perturbation, central- 
body switching, etc. The two-body version does not contain any calibration and only two-body equations 
are used. Therefore no high-fidelity trajectory is obtained. 
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A. Optimality of the optimal two-body abort targeting algorithm 

When the vehicle is located outside the pseudostate sphere, the optimal two-body abort algorithm will be 
used (see line 7 in Algorithm 2 and also in Ref. 4 ). In this section, we will compare the Ac associated with the 
abort maneuvers obtained with the two-body version of the algorithm to the ones obtained using a numerical 
optimizer. The index used for the performance evaluation of the targeter is: 

, |Av (ar — Av*| 

relative error = 100 — 

Av* 

where Av tar is the solution obtained using the two-body targeting algorithm and Av* is the solution 
obtained by the numerical optimizer. 

The numerically optimized solution will be computed taking into account high-order gravity models, 
perturbations, etc. The gravity model used in all the examples consists of an 8x8 Earth and includes 
the Moon and the Sun as point masses. The Copernicus 1,3 trajectory optimization tool has been used to 
numerically optimize the initial guess provided by the targeter. For a given nominal trajectory, different 
impulsive abort maneuvers are calculated for different Abort Times (AT) and Return Times (RT). AT is the 
time elapsed from the Trans-Lunar- Injection maneuver (TLI) to the abort maneuver time. RT is the time 
elapsed from the abort maneuver time to EL As for the El constraints, two cases are considered: (a) equality 
constraints in the El radius and flight-path angle and (b) same as in case (a) plus inequality constraints in the 
El azimuth. In Figure 4, we can see that the relative error between the actual Av and the one calculated by 
the targeting algorithm is always less than 1% for the cases examined. Whereas the relative error increases 
when RT increases (this is due to the fact that the trajectory is affected by third-body perturbations as RT 
increases), the relative error due to different ATs is not as relevant as the one produced by changing RT. This 
is due to the fact that the Av magnitude decreases when AT increases (the vehicle is getting further from 
the Earth) and therefore the relative error remain the same, although the absolute error increases. In Figure 
4 we can see the magnitude of the initial position vector for different ATs. Given the two-body nature of the 
algorithm, it is expected that its performance deteriorates with this distance (i.e. when ATs increases). 

When comparing both cases in Figure 4, we can also see that in case (b) the relative error between 
different ATs is closer than in case (a). Adding the azimuth constraint reduces the number of possible 
solutions and therefore the trajectories obtained by the two-body targeter and the numerical optimizer are 
closer. 

The computational time for this routine is 0.11 ms in average. Given that the routine is mainly analytical, 
adding more constraints to the problem does not substantially increase the computational time. 

B. Optimality of the optimal three-body abort targeting algorithm 

As in the previous section, the optimality of the three-body version of the algorithm has been tested against 
numerically optimized solutions. For a given nominal trajectory to the Moon, several ATs and RTs have been 
evaluated. The process consisted of generating a feasible solution for several ATs and RTs by Algorithm 2 and 
use this solution as the initial guess for the numerically optimized trajectory. The trajectory is propagated 
using two different segments centered at the Moon and at the Earth. The gravity model used for the 
computation of the solutions is: 

• Moon-centered segment: 8x8 Moon, Earth and Sun point masses. 

• Earth-centered segment: 8x8 Earth, Moon and Sun point masses. 

The performance is evaluated based on three parameters: relative and absolute error in Ac and CPU factor: 

relative error 
absolute error 
CPU factor 

where Avtar and Ac* are defined as in the previous section, CPU time optimizer is the time it takes for Coper- 
nicus to find the optimal solution (the overhead associated with Copernicus initialization is not included) 
and C PUtime P seudostate is the time it takes for Algorithm 2 to find the solution. 

8 or 14 


= 100 


, Ac tar — Ac* 


Av* 

Av tar ~ Ac* 

CPU time op timi zer 
CPU tiTfiep Seu( i os t a t e 
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Two-body Abort Targeter Performance 



Two-body Abort Targeter Performance 



(b) 


Figure 4. Comparison between the high-fidelity numerically optimized solution and the two-body solution 
for three different ATs and several RTs. The results are expressed in error Av percentage, (a) El radius 
and flight-path angle equality constraints, (b) El radius, flight-path angle equality constraints and inequality 
constraints in the El azimuth [30°, 60°]. 
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A summary of the test cases evaluated is in Table 1. This table also shows the CPU times associated 
with Algorithm 2 for each of the test cases. We can use the CPU time in this table together with the CPU 
factor presented in all the figures (for example, see Figure 5) to compute the average CPU time associated 
with the numerically optimized solution. 

From the results in Figures 5, 6 and 7, we can conclude that the trajectories generated by Algorithm 2 
are simply feasible trajectories. Although, the calculation of the primary conic is done in an optimal way, 
this does not guarantee the optimality of the solution. The maximum deviation in Av found in the examples 
is about 20% for flyby cases and 6% for direct return cases. But the time it takes for the optimizer to find 
the solution is several times the time it takes for Algorithm 2 (see Figures 5, 6 and 7 to examine specific 
cases). 

It is interesting to note that when the number of constraints in the problem increases (reducing the 
number of possible trajectories) the solutions obtained by Algorithm 2 are closer to the ones obtained by the 
optimizer. A similar effect was shown in the previous case for the two-body version of the abort targeter. As 
an example, when the El Azimuth is constrained as in cases 5 and 6 in Table 1, the maximum relative error 
is reduced to less than 2% (see Figure 7). In Table 1, we can also see that adding constraints to Algorithm 
2 does not affect substantially the CPU time. This is due to the fact that all the El constraints are handled 
in the computation of the primary conic and this is done analytically in a process with no iterations. 


Table 1. 3-Body Targeter (Algorithm 2) Computational Time Performance 


Case 

ALTD E i 

7 EI 

Az E i 

Flyby 


AT=65h 

AT=75h 

AT=85h 

1 

/ 

/ 

X 

X 

mean 

0.34 s 

0.19 s 

0.18s 


std 

0.27 s 

0.05 s 

0.03 s 

2 

/ 

/ 

X 

/ 

mean 

0.23 s 

0.15 s 

0.16 s 


std 

0.03 s 

0.01 s 

0.01 s 

3 

/ 

/ 

[30°, 60°] 

X 

mean 

0.20 s 

0.19 s 

0.18 s 


std 

0.07 s 

0.07 s 

0.06 s 

4 

/ 

/ 

[30°, 60°] 

/ 

mean 

0.31 s 

0.17 s 

0.17 s 


std 

0.03 s 

0.015 s 

0.01 s 

5 

/ 

/ 

[44°, 46° ] 

X 

mean 

0.2 s 

0.19 s 

0.19 s 


std 

0.07 s 

0.06 s 

0.05 s 

6 

/ 

/ 

[44°, 46° ] 

/ 

mean 

0.29 s 

0.18 s 

0.18 s 


std 

0.03 s 

0.02 s 

0.01 s 


V. Conclusions 

Generating an optimal abort return trajectory for Earth to Moon transfers can be computationally 
expensive when considering realistic El constraints and high-order gravity models at the Earth and Moon 
together with third-body perturbations. The use of an initial guess generator can reduce substantially 
the computational time. In this paper, a two-body and three-body abort targeting algorithms have been 
presented. They can both be used as initial guess generators for numerical optimizers. We have shown 
that their computational time is several orders of magnitude less than the ones associated with numerical 
optimizers. We have also discussed the optimality of the solutions found by both targeters. We have 
concluded that whereas the two-body version of the targeter produces solutions close to the optimal (even 
for cases not close to the Earth), the three-body version of the algorithm generates only good feasible 
solutions. Depending on how constrained the problem is, those solutions can be very close to the optimal 
too (see Figure 7). The main advantages of the proposed algorithms are: 

• They can produce feasible solutions in a reasonable amount of time: 0.11 ms for the two-body version 
and 0.5s for the three-body version. It is worth mentioning that the three-body version incorporates 
calibration, which increases the computational time but at the same time allows the generation of 
high-fidelity feasible solutions (including realistic El constraints and gravity fields). 
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Three-body Abort Targeter Performance for non-Flyby Returns 



Three-body Abort Targeter Performance for Flyby Returns 



(b) 


Figure 5. Performance evaluation of the Three-body targeter for Cases 1 (top) and 2 (bottom) in Table 1. 
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Three-body Abort Targeter Performance for non-Flyby Returns 




Three-body Abort Targeter Performance for Flyby Returns 



(b) 


Figure 6. 


Performance evaluation of the Three-body targeter for Cases 1 (top) and 2 (bottom) in Table 1 
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Three-body Abort Targeter Performance for non-Flyby Returns 



(a) 


Three-body Abort Targeter Performance for Flyby Returns 




(b) 


Figure 7. Performance evaluation of the Three-body targeter for Cases 1 (top) and 2 (bottom) in Table 1 
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• The type of solution (flyby or direct) can be specified. Whereas numerical optimizers might pick the 
optimal one of the two approaches, these algorithms allow the user to evaluate a specific one. 

• They can be easily modified to solve some similar one- impulse return trajectory problems. Although 
not shown in this paper, Algorithm 2 has been also used to calculate one-impulse trans-Earth-injection, 
partial-TLI and partial-LOI abort return trajectories with similar performance to the one described 
above. Hardly no modifications were needed to handle these additional problems. 
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